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Electric power-systems are one of the most important critical infrastructures. In recent years, they have been 
exposed to extreme stress due to the increasing demand, the introduction of distributed renewable energy 
sources, and the development of extensive interconnections. We investigate the phenomenon of abrupt 
breakdown of an electric power-system under two scenarios: load growth (mimicking the ever-increasing 
customer demand) and power fluctuations (mimicking the effects of renewable sources). Our results on real, 
realistic and synthetic networks indicate that increasing the system size causes breakdowns to become more 
abrupt; in fact, mapping the system to a solvable statistical-physics model indicates the occurrence of a first 
order transition in the large size limit. Such an enhancement for the systemic risk failures (black-outs) with 
increasing network size is an effect that should be considered in the current projects aiming to integrate 
national power-grids into "super-grids". 

Our society depends on critical infrastructures which deploy us essential goods and services 1 . Most of such 
infrastructures have a network structure. In the last decade, stemming from the availability of large data 
and based on the statistical physicist perspective of the graph theory, a new paradigm to describe large 
networks has blossomed: the field of Complex Networks 2 ' 3 . The Complex Network approach is a powerful and 
unifying tool that enables to treat on the same footing widely different networked systems, ranging from biology 
to sociology to power grids to Internet and the World Wide Web. In fact, such an approach is not tailored to 
reproduce the detailed behaviour of a system, but to highlight the possible occurrence of collective (systemic). 
Hence, the statistical tools of Complex Network Analysis may be of great use to understand salient properties of 
complex human engineered systems that go beyond the single component behaviour. 

Among critical infrastructures, the electrical power grid is perhaps the most crucial one as many other facilities 
like telecommunications, banking systems, oil and gas pumping, and even water depend on the electric power 
system (EPS) 4 . EPSs have a long history and since their first appearance in 1881 at Godalming in England they 
have evolved into one of the most well engineered and robust network infrastructures; nowadays, higher stan- 
dards for the reliability of interconnected EPS are being developed both in Europe and the US by councils and 
associations of EPS operators 5,6 . Nevertheless, not only power outages do occur, but also large outages are more 
likely than what would be naively expected; in fact, the analysis of historical data reveals that their occurrence has 
fat tails 7,8 , implying a significant risk of system- wide failures. Given the disruption and economic damages caused 
by major outages, understanding the nature of such occurrences is a major problem to be addressed. 

An important general question is whether EPSs are subject to emergent behaviour or not. In fact, EPSs are 
aggregations of large number of simple units and can be analysed as a complex network 9 . It therefore makes sense 
to ask if EPSs, as a whole, exhibit additional complexity beyond what is dictated by the simple sum of its parts. To 
this aim, we investigate if an abrupt breakdown transition could emerge in a simple yet realistic complex network 
model of power grids focused on cascade failures. 

In the context of power systems, a cascading outage is a sequence of failures and automatic disconnections 
consequent to an initiating event; a system-wide outage is also called "black-out". The rapid succession of 
automatic reactions in an EPS happens in a time-scale that is typically too short to stop the process by human 
intervention. Reactions following an initiating event or events include sequential tripping (disconnection) of 
transmission lines and generators. Initiating events can be due to natural causes (like a line sagging into vegeta- 
tion, or high wind or lightning shorting a line) but also to human actions (or inaction) or due to imbalances 
between load and generation. 

While no two cascading outages are the same 10 , we will study a class of possible outages and analyse their 
characteristics. In particular, we will consider the fragility of EPSs with respect to outages due to cascades of line 
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overloads causing lines to trip. To this end, we will put under stress 
both model and realistic EPSs to understand the nature of systemic 
outages. The nature of the stress will be twofold: first, we will consider 
the case in which an increasing demand on a fixed infrastructure 
leads to line overloads and subsequent outages. This would corre- 
spond to the case of EPS that are operated to the limit of their 
capacities in order to maximize profits. Second, we will consider 
the important case of fluctuations in demands and generation; this 
is a particularly relevant case since the steady penetrations of renew- 
able sources is introducing in the grids new erratic sources whose 
effects and consequences on existing power grids have not yet been 
fully understood. 

Results 

Model. A particular source of stress to EPS comes to the fact that 
adjustments in power generation are not real-time but follow fixed 
time schedules; for example, in Europe the production is fixed in 
advance the day before and periodic adjustments happen every 15 
minutes; therefore, the reaction time (apart from automatic controls/ 
tripping) is generally much higher than the time of propagation of 
electrical perturbation in the system. 

The tripping of lines and generators above their operating limits 
induced by automatic protective equipments is common to all kinds 
of cascading outages; while this process is intended to protect costly 
equipments from damage, it can potentially widen cascade failures 10 . 

Most of the cascade model for power grids are purely topological 
models based on the local redistribution of power loads upon fail- 
ure 11 " 14 and disregard the long range nature of electricity. On one 
hand, a clear signature of the non-locality of power outages can be 
found in real data: this is, for example, the case of the tripping 
sequence of the Western Interconnection WSCC system distur- 
bances in July 2-3, 1996 15 , where the occurrence of subsequent fail- 
ures in far away lines can be observed. On the other hand, it is 
possible to have simple yet realistic models respecting Kirchoff laws, 
like the one introduced by Carreras et al 16 . In this paper, we will 
consider a very similar model introduced by Pahwa et al 17 . 

In Pahwa's model, the power grid is modelled as a complex net- 
work where nodes represent buses and links represent electrical 
branches; the initial distribution of loads and sources represents 
the stress imposed to the grid. The initial power flows on lines are 
calculated using the DC power flow model (see Materials). If the load 
on a line goes beyond its capacity, the line trips (disconnects) and 
power flows are recalculated on the new topology (i.e. the grid 
MINUS the tripped lines). Such procedure is repeated until conver- 
gence (fig. 1). We will refer to such a model as the Overload Cascade 
Model (OCM). 

Our model does not account for the whole complexity of an EPS; 
in particular, it disregards both transients and the dynamics of phase 
angles. Nevertheless, it allows to sort out the role of a class of events 
always present in any black-out, i.e. line overloads. Another import- 
ant property of the OCM is that, due to the long-range nature of the 
interaction, it is amenable to analytic approximations that lead to 
predict the universal behaviour of the system (see Methods). 

To analyze cascades, we will use both IEEE model grids 1819 and the 
Polish Grid 19 . Notice that IEEE grids are always employed as a real- 
istic test-bed and are often used as sub-units in the design of real 
power grids 18 . 

Load growth. A very effective way to understand the breakdown of a 
system is to stress it up to the rupture point. To such a scope, we 
consider a growth in the power demand while keeping the network 
fixed. Such mechanism is not so far from reality, since in recent years 
the economic competition and deregulation has led the power 
systems to be operating fairly close to their limits: as an example, 
the number of daily transactions in the US has grown by roughly 
400% since 1998 20 . We model such growth of the demand as a 
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Figure 1 | Algorithm for the Overload Cascade Model (OCM). According 
to the initial distribution of loads and sources (representing the stress 
imposed to the power grid), initial power flows are calculated. If the load 
on a line goes beyond its capacity, the line trips (disconnects) and power 
flows are recalculated on the newtopology (i.e. the grid MINUS the tripped 
lines). Such procedure is repeated until convergence. Notice that, while in 
this paper we employ the DC power flow model 31 (see Materials), the OCM 
algorithm is independent from the power-flow model applied to compute 
line loads. 

simultaneous increase of all the loads by a factor a and record the 
fraction /of branches that fail at the end of a cascade. 

Due to the long-range nature of electric interaction, mean-field 
arguments (see Methods) predict the occurrence of a first order 
transition in the system. Fig. 2 shows the results of the increase of 
load in four IEEE model grids. Indeed, we observe that a steep break- 
down of the grids occurs; unfortunately, data sets for only few real 
(istic) networks are available. 

To further investigate the possibility of a first order transition, we 
generate model networks of increasing sizes using the scheme of 
Pahwa et al 21 (see Materials). We observe that each realization of a 
network, due to the randomness, shows a sharp jump in the final 
fraction / of failures at a different value a c of the load increase. 
Therefore, we calculate/by averaging over different networks at fixed 
a — a c . The result for synthetic networks (fig. 3) indeed confirms that 
the transition is first-order. 

Load fluctuations. Another important source of stress in power 
grids are flow fluctuations 22 ; as an example, customer behaviour 
introduces fluctuations that present both cyclic and stochastic 
components. The introduction of renewable energy sources has 
introduced further erratic fluctuations whose size can even invert 
power flow directions. Such a situation can be troublesome in grids 
that have been engineered for a top-down distribution of the power 
and where the flow directions are assumed fixed. In particular, 
fluctuations beyond the operator's predictions can lead to condi- 
tions where power flows exceed the line ratings significantly; in 
general, distributed generation could compromise the stability of 
the grid 23 . 

We parametrize the size of fluctuation by allowing the initial loads 
to fluctuate by a fraction a and calculate the final fraction of cascaded 
links by averaging over load fluctuations. Fig. 4 shows the results of 
the increase of the fluctuations' size a in four IEEE model grids. 
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Figure 2 | Effects of load growth for IEEE grids. Plot of the fraction /of tripped lines versus the relative load increase a mimicking the growth in power 
demand. The results show an abrupt jump in /that is consistent with the mean-field predictions for the model (see Methods). 
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Figure 3 | Effects of load growth in synthetic networks. Plot of the fraction /of tripped lines versus the relative load increase a mimicking the growth in 
power demand. Results are averaged over several ( 100) realizations of the networks for each size. Since each generated network experiences such a jump at 
a different critical value a c , the resulting /s are averaged as a function of a — a c . Error bars are of the order of the size of the symbols. Our results 
show that also in synthetic networks /has an abrupt jump consistent with the mean-field predictions for the model (see Methods). 
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Figure 4 | Effects of load fluctuations for IEEE grids. Plot of the fraction /of tripped lines versus the relative strength of fluctuations cr. Comparing with 
the case of increasing load of fig. 1, we see that randomness smoothens the transition; yet, the transition becomes steeper with size, indicating the 
possibility of a firts order transition in the limit of large networks. Results are averaged over 100 realizations of the noise. Error bars are of the order of the 
size of the symbols. 



Compared to the load increase case, we observe that the breakdown 
of the grids happens in a smoother fashion but becomes steeper with 
increasing system size, indicating the possibility of a first order trans- 
ition in the limit of large networks. Moreover, we observe that a finite 
fraction of the power lines survives past the transition region, i.e. the 
system fragments in "islands". 

To better explore such an effect, we again analyse the behaviour of 
model networks generated according to the scheme of Pahwa et al 21 
(see Materials) and average the results over both fluctuations and 
different networks of ths same size. The results for synthetic net- 
works (fig. 5) indicate that in the case of fluctuating sources there 
are strong finite size effects smoothening the transition; nevertheless, 
the transition becomes steeper with growing system sizes, indicating 
the possibility of a first order transition. Nevertheless, also in this case 
the breakdown due to fluctuations is less severe since a finite fraction 
of power lines survives and the system is in an "islanded" state. 

Polish grid. We now consider a bigger network of 2746 nodes, the 
Polish grid (available in 19 ). This network is a snapshot of the national 
high-voltage Polish power grid obtained from the data collected and 
used by Polish transmission system operators. Such a network is 
realistic and, thanks to its size, helps us to further confirm the 
increasing abruptness of the breakdowns with the size. The upper 
panel of fig. 6 shows that - as expected - by increasing all the loads in 



such a grid, the breakdown is abrupt as in a first order transition. The 
situation is more complicated in the case of random loads. The lower 
panel of fig. 6 shows the un- averaged results for different realization 
of the noise (fluctuations). For such a big network we find that for a 
given realization of the noise the system is either in a safe state (the 
fraction of tripped links is ~0) or in a systemic failure state (the 
fraction of tripped links is ~1). At difference with the case of 
uniform load increase, the transition does not happens at a given 
(7, but can happen in a whole range of values. In particular, for the 
Polish grid we find that the grid stays essentially intact when the loads 
are allowed to fluctuate less than —6% and that the system comes in a 
black-out state when loads are allowed to fluctuate more than —12% 
of their nominal values. For intermediate values, the system can 
either be in a safe state or in a black-out state depending on the 
realization of the disorder. Notice that for such a large network the 
breakdown is systemic (the fraction of failed links /is either ~0 or 
~1) and no "islanding" is present even in the fluctuating case. 

Discussion 

Our analysis of real, realistic and model systems indicates that black- 
outs due purely to line overloads are intrinsically a first order phe- 
nomenon as predicted by a simple mean-field modelling of the 
complex network (see Materials). In particular, stressing the system 
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Figure 5 | Effects of fluctuations in synthetic networks. Plot of the fraction /of tripped lines versus the relative strength of fluctuations cr. Comparing with 
the case of increasing load in synthetic of fig. 2, we see that randomness smoothens the transition; yet, the transition becomes steeper with size, indicating 
the possibility of a firts order transition in the limit of large networks. Results are averaged over 10 networks and over 10 realizations of the noise 
for each network. Error bars are of the order of the size of the symbols. 

increasing the size of the system enhances its vulnerability to sys- 
temic cascades. 

Methods 

Power Flow models. Power flow models aim to obtain complete information on 
voltage angles and magnitudes at each bus of a power system at given loads and 
generation. For an AC (alternate current) system, the power balance equations can be 
written for each bus (nodes of the network); real and reactive power flow on each 
branch (edges of the network) and the generator reactive power output can be 
analytically determined, but due to the non-linear character of the system numerical 
methods are employed to obtain a solution 30 . A bus of the system is either classified as 
Load Bus if there are no generators connected or as a Generator Bus if one or more 
generators are connected. It is assumed that the real power P D and the reactive power 
Q D at each Load Bus are given, while for Generator Buses the real generated power P G 
and the voltage magnitude | V| are given. A particular Generator Bus, called the Slack 
Bus, is assumed as a reference and its voltage magnitude | V| and voltage phase @ are 
fixed. The branches of the electrical system are described by the bus admittance 
matrix Y. 

DC power flow. A simplification of the AC power flow equations is obtained by 
linearising the equations by requiring that bus voltages V f are fixed and that phase 
differences 9^ = © f — @ ; along the branches are small. The resulting linear system of 
equations constitutes the DC power-flow model 31 

P = £@ (1) 



by uniformly increasing the power dispatched in the grid causes a 
systemic failure of the lines. 

On the other hand, when the breakdown is fluctuation- driven, 
such a transition is smoothed at finite sizes: this means that bigger 
systems are more fragile, since their breakdown happens in a more 
abrupt fashion. On the other hand, in the case of fluctuation- driven 
breakdowns, small system are seen to survive in a fragmented state 
where separate "islands" survive: such an effect is not predictable by 
our mean-field approach and is possibly due to the planar nature of 
real and realistic grids. 

Hence, our results indicate that the realization of the foreseen 
super-grid projects 24 ' 25 could introduce an intrinsic source of sys- 
temic risk: in fact, since joining together more networks increases 
the size of the global grid, it could also enhance its fragility. 

By considering static power flow models, we contribute only on a 
single aspect of power grid vulnerability, i.e. the "static security 
assessment" 26 . Hence, we disregard other important phenomena like 
line tripping due to frequency imbalances generated by the mismatch 
between demand and generation 27 . Recent studies hint that tuning 
the parameters of a grid using a complex network approach could 
stabilize its dynamics mitigating the impact of such transients 28 ; in 
general, changing the topology of a network could improve its res- 
ponse to attacks and failures 29 . 

In conclusion, while we have not meant to fully address the whole 
phenomenology of blackouts and outages in power grids, we have 
isolated one of the key factor that concur to such phenomena, show- 
ing by simulations and the means of a simple mean-field model that 
overload failures are intrinsically a first order phenomenon where 



where P f is the total power (generation minus load) at the i-th bus, C = K — Y is the 
Laplacian matrix with K diagonal degree matrix K u = ^\ Yy; the sign of P z deter- 
mines if a node is a generator, a load or even a transit node (P, = 0). Notice that we are 
neglecting phase shifts of the transformers that would add an extra term in eq.(l): 
P = £0 + P^. 

For the purposes of EPSs operators, the DC approximation matches fairly enough 
with the full AC solution 32 ; in fact, DC power flow is on average wrong by a few 
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Figure 6 | Cascade results for the 2746 bus Polish grid. This network is a 
snapshot of the national high-voltage Polish power grid obtained from the 
data collected and used by Polish transmission system operators. In the 
upper panel we show that the fraction /of tripped lines versus the relative 
load increase a shows and abrupt jump typical of first order transitions; 
moreover, breakdowns in such a network are systemic, i.e either all the 
lines have failed (f~ 1) or the grid is safe (f~ 0). In the lower panel, by 
varying cr, we show the effects of the fluctuations in the loads. Data point 
corresponds to 10 different realizations of the noise. At difference with the 
case of uniform load increase, the transition does not happens at a given cr, 
but can happen in a whole range of values. In particular, we find that the 
grid stays essentially intact (f~ 0) when the loads fluctuate less than —6% 
and that the system suffers a black-out (/— 1) when loads fluctuate more 
than —12% of their nominal values. For intermediate values (shaded 
region), the system can either be in a safe state or in a systemic black-out 
state depending on the realization of the disorder. 

percent 33 with respect to the more computationally intensive AC power flow. Another 
important characteristic is that the probability distribution of the line loads looks 
statistically equivalent for the AC as well as the DC approximation 34 . Neglecting line 
resistances is a justified assumption because they are very small as compared to their 
inductive reactances 31 . 

Synthetic power grids. The topology of a network determines the arrangement of 
nodes in the network and how they are connected to each other. A small change in the 
network topology, such as removal or addition of nodes and/or links, may lead to 
changes in the network properties. These changes, in turn, may affect the robustness 
of the network. With the increasing penetration of distributed sources, the network 
structure would play a very important role in deciding the threshold beyond which 
any increase in distributed generation could cause a cascade failure in the system. 

There is limited real data available for power grids. In order to have more networks 
and also to study the effect of different topologies on the robustness of the power grid 
due to increasing amount of distributed generation, it is necessary to have model 
topologies, i.e a network generation algorithm to generate synthetic power grid 
networks. 

The Erdos-Renyi model 35 and the network generation algorithm of Molloy and 
Reed 36 are the simplest network models that include stochasticity as an essential 
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Figure 7 | Synthetic networks of different sizes. Our model networks are 
generated using a variation of the Generalized Random Graph Model that 
imposes constraints on the statistics such that the average degree, the 
maximum degree of the network under construction, the fraction of 
generation nodes, the fraction of load nodes and the generation/load 
statistics are consistent with the ones found in the empirical IEEE test 
networks. 

element in the construction of the network. In these two models the connection of the 
nodes is completely random and lacks knowledge of the principles that guide the 
creation of links between nodes. Thus, links in these models are created randomly 
with a given connection probability p. Another model for creation of random graphs 
is the Generalized Random Graph model 37 " 39 which still keeps the assignment of links 
random but specifies a predefined degree sequence. Thus, the Generalized model can 
produce graphs with degree distributions which are not necessarily Poisson. The 
Generalized Random Graph model was first proposed in 1978 in 40 and is also known 
as the Configuration model. 

We notice that, although some steps in the directions of generating realistic power 
grid topologies have been done 41,42 , there is not yet a consensus on the networks to 
employ as model EPSs; for example, in their studies of inter-dependent power grids, 
Brummit et al. use generalized random regular graphs to generate their model 
networks 43 . 

Network generation. Our model networks are generated using a variation of the 
Generalized Random Graph Model and imposing certain constraints on the 
statistics 21 . Such constraints come in the form of specifying the desired average 
degree, the desired maximum degree of the network under construction, the 
fraction of generation nodes, the fraction of load nodes and the generation/load 
statistics. For each size N, we generated statistically independent configurations 
with the same number of links. This was maintained by considering a constant 
average degree for each of the configurations. The average node degree is given by 
<k> = 2 M/N where M is the number of links and N is the number of nodes. As 
an example, for a network of size N = 50, if the average node degree, <k> is 
considered as 2.8, the number of links M is 70. The impedances, the initial loads 
and the generators are generated by selecting samples from distributions derived 
from the impedance, load and generation data of the available real power grid 
networks 18 . 

In addition, since the generated networks should be as close to realistic networks as 
possible, we also match the maximum node degree of the generated networks with the 
corresponding test networks. Extra care has been taken to ensure that we obtain a 
connected network with no self loops. In fig. 7 we show the topologies of some of the 
synthetic networks. 

Simulations. We stress real and synthetic networks according to the two protocols of 
increasing loads or fluctuations. In the first case, loads are uniformly increased 



L i (a) = (l + a)I? 



(2) 



where L° are the initial values of the loads and a parametrises the strength of the stress. 
In the case of synthetic networks, results are averaged over 100 samples for each 
system size. 

In the second case, fluctuations are mimicked by using a random variable £ f uni- 
formly distributed in [ — 1, 1]: 



SCIENTIFIC REPORTS | 4 : 3694 | DOI: 1 0.1 038/srep03694 



6 



If 



0.001 



H 1 1 I I I I I | 1 1 1 I I I I I | E 

14 bus 




— y=-0.01 



i.01 



_i I I I I I 



0.1 



If 



0.1 



- 0.01 - 



0.001 



~i 1 — i i i i i 1 1 1 1 — i i i i i 1 1 

30 bus 




— y=-0.7 



00 



i.01 



_l I I I I I 



0.1 



0.1 

P(c) 

0.01 



i — i — i 1 1 1 1 1 1 
50 bus 



01 




0.1 



- 0.01 - 



000 o. 



~l 1 I I I I I I 1 1 1 I I I I I I 

247 bus 




: — Y=-l-8 o ©00 O 



J I I I I I I I I I I I LLXJJJ 



01 



0.1 



Figure 8 | Histograms of the normalised line capacities cHere c = C/C max and C max is the maximum line capacity of the network) for the 14, 30, 50 and 
247 bus IEEE networks. The thick lines indicates fits of the form p{c) oc c~ 7 . Notice that the exponent 7 increases with system size. 



(3) 



i.e. Li{o) e [(1 — o)L b (1 + o)Li]. Here o parametrises the strength of the stress. In 
the case of the IEEE networks and of the Polish grid, averages are obtained over 100 
realizations of the noise. In the case of synthetic networks, averages are obtained over 
10 realizations of the noise for each of 10 separate network samples. 

Mean field analysis. Flows in power grids have complex dynamics even in the DC 
approximation: if a single line gets overloaded or breaks, its power is immediately 
distributed, not to a single different branch, but in the whole system, due to the long- 
range nature of Kirchoff s laws. The redistribution of power after a component failure 
is dependent on the electrical characteristics (such as impedances) of the lines; if the 
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Figure 9 | Mean-field results for increasing loads The fraction of failures f 
grows as a function of the initial load per link until the first order 
transition. For our mean-field model no islanding is possible since the 
fraction of failed nodes after the transition is/= 1 (all the links have failed). 
Left panel: For an initial distribution p(C) oc const, the approach to the 
first order jump is smoother. Right panel: For an initial distribution p( C) 
oc C~ 2 (typical of larger networks), the first order jump happens with no 
apparent precursors. 



redistribution of power leads to subsequent overloadings, the consequence could be a 
cascade of overloading failures. In general, after an initial failure, some of the lines 
could get overloaded and fail: this represents the first stage of cascade. First stage 
could possibly lead to further overloading and collapses, constituting the second stage 
and so on. In this way, the system goes through multiple stages of cascades until it 
finally stabilizes and there are no more failures. 

We indicate the final stable configuration of the OCM by {y}, where 



yr 



0 if branch is broken 
Yn otherwise 



(4) 



and \yij6jj\ < Qp where 6y is the phase angle difference among nodes i,j and Cy is the 
threshold flow (the capacity) beyond which a branch becomes overloaded and fails. 

Following Zapperi et al. 44 , the stable configuration of the OCM can be obtained 
from a double minimization of an objective function 



(5) 



where Cy = Cy/Yy and)/y <E {0, Yyj. Here H is minimized both with respect to the 
global variables 9% at fixed (notice that DC power flow implicitly minimizes the 
"dissipations" yffij) and with respect to the local variables y t j at fixed 0 y - (notice that 
this term corresponds to the breakdown of overloaded links). In such a form, the 
OCM can be mapped in the model for the breakdown of a disordered media of 
reference 44 , indicating that the cascade transition in power grids is a first order 
transition, i.e. consists of an abrupt failure of the system. 

The first-order nature of such a transition stems from the interplay of the long 
range nature of Kirchoff s laws and the distribution of failure thresholds 45 . Long- 
range interactions allow to apply mean field analysis to the model; in particular, it is 
possible to assume that when a line is tripped, its load could re-distribute equally 
among all other lines. Subsequently, the lines above their threshold would trip again, 
their loads would be re-distributed and so on up to convergence: such a model 
corresponds to the so called democratic fiber-bundle model (DFBM) 46,47 . Depending 
from the distributions of threshold, the DFBM can have at the same time an abrupt 
breakdown of the system and a power-law distribution of cascade sizes 45 ; such 
behaviour resembles the hybrid transition for multiplex networks 48 . 

The DFBM model shares the same mean-field spirit of the CASCADE model for 
black-outs 49 ; notice that while DFBM assumes a starting situation in which loads are 
the same and thresholds are randomly distributed, CASCADE assumes that loads are 
randomly distributed and thresholds are the same. While the CASCADE model 
predicts a self- organised second order critical point, our OCM predicts a first order 
transition; a recent analysis on black-outs' sizes in several countries indicates that 
causal factors other than self-organization or a critical state might be significantly 
ruling the system dynamics 8 . 



SCIENTIFIC REPORTS | 4 : 3694 | DOI: 1 0.1 038/srep03694 



7 



{C[link} is the capacity of the link} 
{F[link\ is the initial flow on the link} 
{ M is the initial number of links} 
Msafe = M 
for link in 1 ... M do 

broken[link] ^— false 
repeat 

Stress <— 0 

for all broken[link] = false do 
if F[link] > C[link] then 
Stress ^— Stress + F[link] 
broken[link] <— trite 
Msafe <- Msafe - 1 
for all broken[i] = false do 

F[link] <- F[link] + Stress/Msafe 
until Msafe = 0 or Stress = 0 

return / <- (M - Msafe)/M 

Figure 10 | Algorithm for mean-field with fluctuations in the initial 
loads. 



Finally, we recall that, by identifying phase differences with voltage drops, the DC 
power flow equations are equivalent to the equations for DC currents on a resistor 
network; therefore, our model can also be mapped in the so-called random-fuse 
model 50 . 

In the following, we analise the case of increasing loads and of random fluctuations 
in a simple DFBM-like model for the power grid, showing that mean-field predicts a 
first order transition both in the case of increasing loads and increasing fluctuations. 

DFBM with increasing loads. We consider M lines where each link i is characterised by 
its capacity Q (maximum amount of flow) and assume that capacities are charac- 
terised by a probability distribution p(C). When the load of a link surpasses its 
capacity, the link breaks and its load gets redistributed "democratically" among the 
remaining surviving links. 

We first consider a uniform stress situation, where the system is subject to a load L 
and all the links are initially subject to the same load / = L/M. Let's consider a cascade 

/L/M 
p(C)dC immediately 
o 

fails since their thresholds are less than the load / they sustain. Hence, after the first 
stage of a cascade, there are M 1 = ( 1 — f)M surviving links and the new load per link 
is l l = L/M 1 . The cascade's stages follow analogously; we can thus write the mean field 
equations for the (t + l)-th stage of the cascade: 



l t+1 =L/M f ;l°=L/M 

f +1 =p(i t+l yj° = i 

< M t+1 = (l-f t+1 )M;M°=M 



(6) 



where P(x)= p(C)dC is the cumulative distribution function associated to p(C). 
Jo 

Eq. 6 can be simplified in a single equation for the fraction /of failed links 



0.8 



• Simulation 
— Mean Field 




0.3 



0.4 



0.5 



(7) 



Figure 1 1 | Mean-field results for increasing fluctuations. When we 
consider fluctuations A =7 in the initial loads around their average value 7 
we observe that the transition is still abrupt but can happen in a region 
around the mean-field predictions (corresponding to A = 0). Full circles: 
results of simulation for 1000 links systems where, for each value of /, 10 
different realizations of the disorder in the initial loads are considered. 
Thick line: results for the mean field case where no fluctuations are present 
and all the initial loads are equal to 7. The shaded area highlights the region 
(roughly centred around the mean field transition) where the first order 
transition can occur. 



Hence, the fixpoint/* = F(J*) of eq. 7 represents the fraction of links survived to the 
failure cascade. The behaviour of/* depends on the functional form of p(C) 45 and is 
known to present a first order transition for a wide family of curves. Thus, it is 
important to understand the functional form of the distribution for the thresholds C 
in the case of realistic power networks. In fig. 8 we show the histograms of the 
normalised line capacities c= CI C max where C max is the maximum line capacity of the 
network. For all the networks we find that such histograms can be approximated by 
p's of the form 

p(c)~c~^ for c min <c<l 

where c min ~ 0.02 for all the networks and the exponent 7 increases according to the 
size. 

To highlight the effects of the shape of p(c), we plot in fig. 9 the results for the MF 
model (the fixpoint of equation 7) for the extreme cases of a flat distribution (7 = 0) 
and of a sharp decrease in the distribution (7 = 2). We find in both cases that the 
transition is first order. Moreover, the height of the jump increases for bigger 7> 
indicating that the size behaviour of cascade in power networks could be mainly due 
to the change in shape of the capacities' distribution with size. 

DFBM with increasing fluctuations in loads. To account for the fluctuations in the 
initial loads, we can generalise the DFBM model to the case where the initial loads on 
the links are not uniform. Hence, we consider a model in which M links with capa- 
cities Q extracted with probability p(C) have initial random loads /? and where the 
loads of failed links get redistributed equally among the remaining links. The fraction 
/of unbroken links can be solved numerically by iterating the algorithm of fig. 10. 

The average initial load is 1=M~ 1 /,•; for simplicity, we will consider the case 
where the initial loads f? are uniformly distributed among 0 and 27, i.e. between 
[7— A, 7+ A] with A =7. By solving numerically such system at finite M for different 
values of the noise A = 7 and several realizations of the initial loads, we find that the 
transition does not happen anymore at a critical value of 7 but can occurr in a region 
around the mean-field prediction. As an example, we show in fig. 1 1 the results for 10 
realizations of the initial loads in a system of M = 1000 links. 
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